ANÁLISIS DE MOVILIDAD EN BOGOTÁ

knitr::opts_chunk$set(warning = FALSE, message = FALSE)
options(warn = -1)
setwd("C:/Users/Angelica/Documents/Movilidad/Data")
library(readxl)
Hogares <- read_excel("a. Modulo hogares.xlsx")
#Vehiculos <- read_excel("b. Modulo vehiculos.xlsx")
#Personas <- read_excel("c. Modulo personas.xlsx")
#Viajes <- read_excel("d. Modulo viajes.xlsx")
#Etapas <- read_excel("e. Modulo etapas.xlsx")

ESTADÍSTICAS DESCRIPTIVAS DE HOGARES

UTAM: Unidades Territoriales de Análisis de Movilidad

Cargamos el archivo shapefile

library(sf)
library(ggplot2)
setwd("C:/Users/Angelica/Documents/Movilidad/UTAM")
UTAM_shape <- st_read("UTAM2023.shp")

## Reading layer `UTAM2023' from data source 
##   `C:\Users\Angelica\Documents\Movilidad\UTAM\UTAM2023.shp' using driver `ESRI Shapefile'
## Simple feature collection with 142 features and 6 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -74.44978 ymin: 3.73103 xmax: -73.86269 ymax: 5.041663
## Geodetic CRS:  MAGNA-SIRGAS

Calculamos el conteo de hogares

El codigo UTAM del shapefile es diferente al de hogares dado que, por ejemplo, en hogares la UTAM1 se escribe UTAM001. Por tanto es necesario modificar la variable “cod_utam_hg”

library(dplyr)
Hogares <- Hogares %>%
  mutate(cod_utam = gsub("^(UTAM|UPR)0*", "\\1", cod_utam_hg))

num_hogares <- as.data.frame(table(Hogares$cod_utam))
colnames(num_hogares) <- c("UTAM", "Hogares")  

Luego, lo unimos con el shapefile y generamos el mapa

shape_joined <- UTAM_shape %>%
  left_join(num_hogares, by = "UTAM")

ggplot() +
  geom_sf(data = shape_joined, aes(fill = Hogares), color = "black") +
  scale_fill_viridis_c(option = "C", na.value = "gray90") +
  theme_minimal() +
  labs(title = "Número de hogares por UTAM", fill = "Hogares")

Si se excluyen las zonas de planeamiento rural (UPR) podemos ver más detalladamente las UTAM centrales

shapefile_utam <- shape_joined %>%
  filter(grepl("^UTAM", UTAM))

ggplot() +
  geom_sf(data = shapefile_utam, aes(fill = Hogares), color = "black")+
  scale_fill_viridis_c(option = "plasma", na.value = "grey90") +
  theme_minimal() +
  labs(title = "Número de hogares por UTAM",
       fill = "Hogares")

Para solucionar las dimensiones del mapa, se puede hacer un gráfico interactivo con plotly

library(stringr)
shapefile_utam <- shape_joined %>%
 mutate(UTAMNombre = if_else(str_trim(str_to_upper(UTAMNombre)) == "N/A", MUNNombre, UTAMNombre))%>%
 mutate(label = paste0("Codigo : ", UTAM, "\nNombre UTAM: ", UTAMNombre, "\nHogares: ", Hogares ))
  
mapa <- ggplot() +
  geom_sf(data = shapefile_utam, aes(fill = Hogares, text = label), color = "black") +
  scale_fill_viridis_c(option = "plasma", na.value = "grey90") +
  theme_minimal() +
  labs(title = "Número de hogares por UTAM", fill = "Hogares")

library(plotly)
ggplotly(mapa, tooltip = "text")

Visto en diagrama de barras para municipios diferentes a Bogotá

library(forcats)
ggplot(subset(Hogares, nom_mpio_hg != "Bogotá D.C."), 
       aes(x = fct_rev(fct_infreq(nom_mpio_hg)))) +
  geom_text(stat = "count", aes(label = ..count..), hjust = -0.2, size = 3) +
  geom_bar(fill = "steelblue") +
  coord_flip() +
  labs(
    title = "Número de hogares por municipio",
    x = "Municipio",
    y = "Cantidad de hogares"
  ) +
  theme_bw()

Variable ingresos

ggplot(Hogares, aes(x = ingre_mes_hg)) +
  geom_bar(fill = "steelblue") +
  geom_text(stat = "count", aes(label = ..count..), vjust = -0.3, size = 3) +
  labs(
    title = "Distribución de ingreso mensual",
    x = "Ingreso mensual (COP)",
    y = "Número de hogares"
  ) +
  theme_bw() +
  theme(axis.text.x = element_text(angle = 90, hjust = 1))

# Contar hogares por municipio e ingreso
tabla_prop <- Hogares %>%
  filter(ingre_mes_hg != "NS/NR") %>%
  count(nom_mpio_hg, ingre_mes_hg) %>%
  group_by(nom_mpio_hg) %>%
  mutate(prop = n / sum(n)) %>%
  ungroup()

# Heatmap con proporciones
ggplot(tabla_prop, aes(x = ingre_mes_hg, y = fct_infreq(nom_mpio_hg), fill = prop)) +
  geom_tile(color = "white") +
  scale_fill_viridis_c(option = "C", labels = scales::percent) +
  labs(
    title = "Distribución de ingreso mensual por municipio",
    x = "Ingreso mensual",
    y = "Municipio",
    fill = "Porcentaje"
  ) +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 90, hjust = 1))

ggplot(subset(Hogares, nom_mpio_hg == "Bogotá D.C."),
       aes(x = estrato_hg, fill = factor(ingre_mes_hg))) +
  geom_bar(position = "fill") +
  scale_y_continuous(labels = scales::percent_format()) +
  scale_fill_manual(values = colores) +  # Usamos scale_fill_manual con los colores generados
  labs(
    title = "Distribución de ingreso mensual por estrato (Bogotá)",
    x = "Estrato",
    y = "",
    fill = "Ingreso mensual"
  ) +
  theme_bw()

Quitando a categoría NS/NR

  hogares_bog <- Hogares %>%
  filter(nom_mpio_hg == "Bogotá D.C.", ingre_mes_hg != "NS/NR")

hogares_bog %>%  
  ggplot(
    aes(x = estrato_hg, fill = factor(ingre_mes_hg))) +
  geom_bar(position = "fill") +
  scale_y_continuous(labels = scales::percent_format()) +
  scale_fill_manual(values = colores) +  # Usamos scale_fill_manual con los colores generados
  labs(
    title = "Distribución de ingreso mensual por estrato (Bogotá)",
    x = "Estrato",
    y = "",
    fill = "Ingreso mensual"
  ) +
  theme_bw()

O visto no en proporcion

ggplot(subset(hogares_bog, nom_mpio_hg == "Bogotá D.C."),
       aes(x = ingre_mes_hg, fill = factor(estrato_hg))) +
  geom_bar(position = "dodge") +  # Posición dodge para que las barras estén separadas
  labs(
    title = "Distribución de ingreso mensual por estrato (Bogotá)",
    x = "Estrato",
    y = "Número de hogares",
    fill = "Ingreso mensual"
  ) +
  theme_bw()

Tamaño de los hogares

ggplot(Hogares, aes(x = perstotal_hg)) +
  geom_histogram(binwidth = 1, fill = "steelblue", color = "white") +
  labs(
    title = "Distribución tamaño del hogar",
    x = "Número de personas en el hogar",
    y = ""
  ) +
  theme_bw()

ESTADÍSTICAS DESCRIPTIVAS DE VEHICULOS

Tipos de vehículo

Se cuenta con una categorización bastante amplia con frecuencias bajas para algunos tipos de vehículo

kable(tabla_vehiculos, caption = "Tipos de vehículo")

Tipos de vehículo

Tipo de Vehículo	Cantidad
Automóvil	6441
Automóvil o camioneta de servicio especial	64
Bicicleta motor eléctrico	55
Bicicleta motor gasolina	9
Bicicleta niños	329
Bicicleta sin motor	4691
Bicitaxi motor eléctrico	9
Bicitaxi motor gasolina	13
Bicitaxi sin motor	190
Camión	83
Campero/Camioneta	901
Moto - carro	134
Motocicleta	4028
Otro	69
Patineta	70
Patineta eléctrica	92
Pick Up/Van	69
Taxi	115
Triciclo - moto	3
Vehículo tracción animal	3
Vehículo tracción humana	24

Por tanto, es necesario realizar una recategorización de la variable. Se propone la siguiente manera

Recategorización de tipos de vehículo

Categoría.nueva	Tipos.originales.incluidos
Automóvil	Automóvil
Motocicleta	Motocicleta
Campero/Camioneta	Campero/Camioneta
Otros motorizados	Taxi, Automóvil o camioneta de servicio especial, Moto - carro, Camión, Pick Up/Van, Triciclo - moto, Bicicleta motor gasolina, Bicitaxi motor gasolina, Bicicleta motor eléctrico, Bicitaxi motor eléctrico, Patineta eléctrica
Otros	Bicicleta sin motor, Bicicleta niños, Bicitaxi sin motor, Patineta, Vehículo tracción humana, Vehículo tracción animal

ggplot(Vehiculos, aes(x = fct_rev(fct_infreq(tipo)))) +
  geom_bar(aes(fill = tipo)) +
  labs(
    title = "Distribución tipos de vehículos",
    x = "",
    y = "Cantidad de vehiculos"
  ) +
  theme_bw() +
  coord_flip() +
  scale_fill_brewer(palette = "Set2") +
  theme(legend.position = "none")

Tipo de combustible

Se propone una recategorización de la variable

Recategorización tipo de combustible

Categoría nueva	Categoria Original
Gasolina	Sólo Gasolina
Diésel	Diésel
Eléctrico	Eléctrico
No aplica	No aplica
Otro	GNV, GNV y gasolina, Híbrido (eléctrico-gasolina / diésel)

ggplot(Vehiculos, aes(x = fct_rev(fct_infreq(combustible)))) +
  geom_bar(aes(fill = combustible)) +
  geom_text(stat = "count", aes(label = ..count..), hjust = -0.2, size = 3) +
  labs(title = "Distribución por tipo de combustible",
       x = "Tipo de combustible", y = "Cantidad de vehiculos")+
  theme_bw() +
  coord_flip() +
  scale_fill_brewer(palette = "Set2") +
  theme(legend.position = "none")

ggplot(Vehiculos, aes(x = (fct_infreq(tipo_placa)))) +
  geom_bar(aes(fill = tipo_placa)) +
  labs(title = "Distribución por tipo de placa",
       x = "Tipo de placa", y = "Cantidad de vehiculos")+ 
theme_bw() +
  scale_fill_brewer(palette = "Set2") +
  theme(legend.position = "none")

UTAM

Grafico excluyendo zonas de planeamiento rural UPR

Vehiculos %>%
  ggplot(aes(x = factor(estrato_hg))) +
  geom_bar(aes(fill = tipo),position = "fill") +  
  scale_y_continuous(labels = percent_format()) +
  labs(
    title = "Distribución del tipo de vehículo por estrato",
    x = "Estrato",
    y = "Proporción",
    fill = "Tipo de vehículo"
  ) +
theme_bw()

library(tidyr)
#Filtrar Vehículos motorizados
vehiculos_filtrados <- Vehiculos %>%
  filter(tipo != "Otros")

#Contar vehículos por hogar
hog_vehiculos <- vehiculos_filtrados %>%
  group_by(cod_hog) %>%
  summarise(num_vehiculos = n())

#Unir con Hogares y reemplazar NA por 0
Hogares_vehiculos <- Hogares %>%
  left_join(hog_vehiculos, by = "cod_hog") %>%
  mutate(num_vehiculos = replace_na(num_vehiculos, 0))

ggplot(Hogares_vehiculos, aes(x = factor(num_vehiculos))) +
  geom_bar(fill = "steelblue") +
  labs(
    title = "Distribución de número de vehículos motorizados en el hogar",
    x = "Cantidad de vehículos",
    y = "Número de hogares"
  ) +
  theme_bw()

Creo que en hogares cuyo ingreso es de 0-0.4M debe haber bicitaxis que generan ese incremento en el numero de vehiculos porque no me parece normal